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Abstract 

We present a complete statistical analysis of quantum optical measurement 
schemes based on photodetection. Statistical distributions of quantum observ- 
ables determined from a finite number of experimental runs are characterized 
with the help of the generating function, which we derive using the exact 
statistical description of raw experimental outcomes. We use the developed 
formalism to point out that the statistical uncertainty results in substantial 
limitations of the determined information on the quantum state: though a 
family of observables characterizing the quantum state can be safely evalu- 
ated from experimental data, its further use to obtain the expectation value 
of some operators generates exploding statistical errors. These issues are 
discussed using the example of phase-insensitive measurements of a single 
light mode. We study reconstruction of the photon number distribution from 
photon counting and random phase homodyne detection. We show that uti- 
lization of the reconstructed distribution to evaluate a simple well-behaved 
observable, namely the parity operator, encounters difficulties due to accumu- 
lation of statistical errors. As the parity operator yields the Wigner function 
at the phase space origin, this example also demonstrates that transformation 
between various experimentally determined representations of the quantum 
state is a quite delicate matter. 
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1 



Typeset using REVTpjX 



I. INTRODUCTION 



Over recent years, the set of tools for measuring quantum statistical properties of optical 
radiation has substantially enlarged. The experimental demonstration of optical homodyne 
tomography has been followed by detailed studies of this beautiful technique and 
diverse novel schemes have been proposed PHM]. The quantum optical "toolbox" for mea- 



suring light contains now experimentally established schemes for reconstructing various rep- 



resentations of its quantum state: the Q function |15[], the Wigner function, and the density 
matrix in the quadrature and the Fock bases A device that is used in most of quantum 

optical schemes to convert the quantum signal to a macroscopic level is the photodetector. 
Thus, photodetection is a basic ingredient of quantum optical measurements. 

The quantum state can be characterized using various representations: quasidistribution 
functions or a density matrix in a specific basis. From a theoretical point of view, all these 
forms are equivalent. Each representation contains complete information on the quantum 
state, and they can be transformed from one to another. Any observable related to the 
measured systems can be evaluated from an arbitrary representation using an appropriate 
expression. 

This simple picture becomes much more complicated when we deal with real experimen- 
tal data rather than analytical formulae. Each quantity determined from a finite number 
of experimental runs is affected by a statistical error. Consequently, the density matrix or 
the quasidistribution function reconstructed from the experimental data is known only with 
some statistical uncertainty. This uncertainty is important when we further use the recon- 
structed information to calculate other observables or to pass to another representation. 
The crucial question is, whether determination of a certain representation with sufficient 
accuracy guarantees that arbitrary observable can be calculated from these data with a rea- 
sonably low statistical error. If this is not the case, the reconstructed information on the 
quantum state of the measured system turns out to be somewhat incomplete. Furthermore, 
transformation between various representations of the quantum state becomes a delicate 
matter. 

These and related problems call for a rigorous statistical analysis of quantum optical 
measurements. The purpose of this paper is to provide a complete statistical description of 
a measurement of quantum observables in optical schemes based on photodetection. This 
approach fully characterizes statistical properties of quantities determined in a realistic mea- 
surement from a finite number of experimental runs. It can be applied either to determi- 
nation of a single quantum observable, or to the reconstruction of the quantum state in a 
specific representation. Within the presented framework we study, using a simple example, 
the completeness of the experimentally reconstructed information on the quantum state. 
We demonstrate the pathological behavior suggested above, when the reconstructed data 
cannot be used to calculate some observables due to rapidly exploding statistical errors. 

An important issue is the statistical methodology applied to retrieve information on the 
measured quantum system from experimental data. Various approaches have been recently 
developed, based on the maximum entropy principle [|16j, the least-squares inversion ||17|| , 
and the maximum likelihood estimation |l8|,|l9j. In this paper we will consider the most 
straightforward and so far the most commonly used strategy, where experimental relative 
frequencies are taken as estimates for quantum mechanical probability distributions. Quan- 
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turn observables are reconstructed from these data via linear transformations. 

Selected aspects of statistical uncertainty resulting from the finite character of the exper- 
imental data sample have been discussed in some particular cases. Systematic and statistical 
errors in homodyne measurements of the density matrix were studied ]5]||, and it was found 
that homodyne detection of quantum observables is accompanied by excess noise compared 
to a direct measurement 0. Analysis of a photon counting scheme for sampling quantum 
phase space showed that compensation for the nonunit detection efficiency is in general not 



possible | [2 Of . The present paper provides a general statistical analysis of quantum optical 
schemes based on photodetection. 

This paper is organized as follows. The starting point of our analysis is the probability 
distribution of obtaining a specific histogram from N runs of the experimental setup. This 
basic quantity determines all statistical properties of quantum observables reconstructed 
from a finite sample of experimental data. We characterize these properties using the gener- 
ating function, for which we derive an exact expression directly from the probability distri- 
bution of the experimental outcomes. These general results are presented in Sec. |I[ Then, 
in Sec. [Til], we use the developed formalism to discuss the reconstruction of the photon 
number distribution of a single light mode, and its subsequent utilization to evaluate the 
parity operator ft. We consider two experimental schemes: direct photon counting using an 
imperfect detector, and homodyne detection with random phase. In both the cases we find 
that the evaluation of the parity operator from the reconstructed photon statistics is a very 
delicate matter. For photon counting of a thermal state, we show that neither the statistical 
mean value of ft nor its variance have to exist when we take into account arbitrarily high 
count numbers. For random phase homodyne detection, the statistical error of the parity 
operator is an interplay of the number of runs N and the specific regularization method used 
for its evaluation. The example of the parity operator illustrates difficulties related to the 
transformations between various experimentally determined representations of the quantum 
state, as the parity operator yields, up to a multiplicative constant, the Wigner function at 
the phase space origin. Finally, Sec. [TV] summarizes the paper. 



II. STATISTICAL ANALYSIS OF EXPERIMENT 

In photodetection measurements, the raw quantity delivered by a single experimental run 
is the number of photoelectrons ejected from the active material of the detectors. The data 
recorded for further processing depends on a specific scheme. It may be just the number 
of counts on a single detector, or a difference of photocounts on a pair of photodetectors, 
which is the case of balanced homodyne detection. It may also be a finite sequence of 
integer numbers, e.g. for double homodyne detection. We will denote in general this data 
by n, keeping in mind all the possibilities. 

The experimental scheme may have some external parameters 6, for example the phase 
of the local oscillator in homodyne detection. The series of measurements are repeated 
for various settings 9{ of these parameters. Thus, what is eventually obtained from the 
experiment, is a set of histograms {/cn}^, telling in how many runs with the settings 9i the 
outcome n has been recorded. We assume that for each setting the same total number of N 
runs has been performed. 
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The theoretical probability p n (@i) of obtaining the outcome n in a run with settings 
$i is given by the expectation value of a positive operator valued measure (POVM) p n (9i) 
acting in the Hilbert space of the measured system. The reconstruction of an observable A 
is possible, if it can be represented as a linear combination of the POVMs for the settings 
used in the experiment: 

^EE^i^Pn^), (1) 

i n 

where a n (8i) are the kernel functions. The above formula allows one to compute the quantum 
expectation value (A) from the probability distributions p n (0i)- 

This theoretical relation has to be applied now to the experimental data. The simplest 
and the most commonly used strategy is to estimate the probability distributions p n (9i) 
by experimental relative frequencies {k n /N)g.. The relative frequencies integrated with the 
appropriate kernel functions yield an estimate for the expectation value of the operator A. 
For simplicity, we will denote this estimate just by A. Thus, the recipe for reconstructing 
the observable A from experimental data is given by the counterpart of Eq. ([J): 

^ = EE^)(|) • (2) 

We will now analyse statistical properties of the observable A evaluated according to 
Eq. (HI) from data collected in a finite number of experimental runs. Our goal is to charac- 
terize the statistical distribution w(A) defining the probability that the experiment yields 
a specific result A. The fundamental object in this analysis is the probability V({k n };6) 
of obtaining a specific histogram {k n } for the settings 9. In order to avoid convergence 
problems, we will restrict the possible values of n to a finite set by introducing a cut-off. 
The probability V({k n };6) is then given by the multinomial distribution: 

where prim in sums and products denotes the cut-off. 

Let us first consider a contribution A, t to the observable A calculated from the histogram 

Bi-. 

4 = I^)(^) fl - (4) 

Its statistical distribution w(Ai] $i) is given by the following sum over all possible histograms 
that can be obtained from N experimental runs: 

w(A t ; di) = E V({k n }; di)5 (A - ^'ani^kn) ■ (5) 

Equivalently, the statistics of Ai can be characterized by the generating function w(X; 9i) for 
the moments, which is the Fourier transform of the distribution w(Af, 6j): 
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w(X;9 l ) = I di i e aii tt)(A i ;^) 



= £ V{{knYA)exp (^E' a n(^)^) • 

{fc n } V " / 

An easy calculation yields the explicit form of the generating function: 



(6) 



w 



(\;6 l )= (l + J2'p n (e i )(e 

\ n 



N 



i\a n (6i)/N 



-V, 



(7) 



The observable A is obtained via summation of the components A{ corresponding to all 
settings of the external parameters Q\. As these components are determined from disjoint 
subsets of the experimental data, they are statistically independent. Consequently, the 
generating function w(X) for the moments of the observable A is given by the product: 



w 



(A) = JdA e iXA w(A) = l[w(\; 8 { 



N 



II l + £P»(fi)(e 



i\a n {9i)/N _ 



1) 



(8) 



This expression contains the complete statistical information on determination of the ob- 
servable A from a finite number of runs of a specific experimental setup. The measuring 
apparatus is included in this expression in the form of a family of POVMs p n {6i)- The 
quantum expectation value of these POVMs over the state of the measured system yields 
the probability distributions p n (0i). Finally, the coefficients a n {6i) are given by the compu- 
tational recipe for reconstructing the observable A from the measured distributions. 

The basic characteristics of statistical properties of the observable A is provided by the 
mean value E(A) and the variance Var(A). These two quantities can be easily found by 
differentiating the logarithm of the generating operator: 



E(A) 
Var(A) 



log w(X) 



l_d_ 
7dA 

1 d 2 



EE' a "(^)Pr, 



(9) 



A=0 i n 



A=0 



1 

N 



i n 



i \ n 



(10) 



The statistical error is scaled with the inverse of the square root of the number of runs N. 
Let us note that the second component in the derived formula for Var(v4) differs from that 
used in the discussions of homodyne tomography in Refs. where it was equal just to 

[E(v4)] 2 . This difference results from different assumptions about the local oscillator phase: 
in Refs. ||,[?| it was considered to be a stochastic variable in order to avoid systematic errors, 
whereas we have assumed that the number of runs is fixed for each selected setting of the 
external parameters. 

The goal of quantum state measurements is to retrieve the maximum amount of informa- 
tion on the quantum state available from the experimental data. Therefore the experimental 
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histograms are usually processed many times in order to reconstruct a family of observables 
characterizing the quantum state. Of course, quantities determined from the same set of 
experimental data are not statistically independent, but in general may exhibit correlations. 
The analysis presented above can be easily extended to evaluation of any number of observ- 
ables from the experimental data. If we restrict our attention to the basic, second-order 
characterization of these correlations, it is sufficient to discuss determination of two observ- 
ables. Let us suppose that in addition to A, another observable B has been calculated from 
the histograms {k n }e i according to the formula: 

^EE'wfe) • (ii) 

in \ iv / 6i 

The generating function w(X, /x) corresponding to the joint probability distribution w(A, B) 
can be found analogously to the calculations presented above. The final result is: 



w(\li) = J dAdB e iXA+i » B w(A,B) 

= II (l + ^Vn(^)(e iAa " (ei)/7V+ ^' l(ei)/JV - 1)) • (12) 
The covariance between the experimentally determined values of A and B is given by: 

Cov( - 4 ' B) = ^ log ' 5(A '' l) 



1 E 

N ~ 



A,/i=0 




(13) 



The covariance can be normalized to the interval [—1, 1] using Var(A) and Var(5), which 
yields the correlation coefficient for the pair of observables A and B: 

Corr(A,B)= Co ^ B ) (14 ) 
/Var(A)Var(5) 



This quantity defines whether the statistical deviations of A and B tend to have the same or 
opposite sign, which corresponds respectively to the positive or negative value of Cott(A, B). 

We have assumed that the histograms k n have been measured for a finite number of 
external parameters settings 9i, which is always the case in an experiment. However, in 
some schemes the measurement of histograms is in principle necessary for all values of a 
continuous parameter. For example, in optical homodyne tomography the full information on 
the quantum state is contained in a family of quadrature distributions for all local oscillator 
phases. Restriction to a finite set of phases introduces a systematic error to the measurement 
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III. PHASE INSENSITIVE DETECTION OF A LIGHT MODE 



We will now apply the general formalism developed in the preceding section to the 
reconstruction of phase-independent properties of a single light mode. The basic advantage 
of this exemplary system is that it will allow us to discuss, in a very transparent way, 
pathologies resulting from the statistical uncertainty. All phase-independent properties of 
a single light mode are fully characterized by its photon number distribution p v . Therefore, 
it is sufficient to apply a phase-insensitive technique to measure the photon statistics of the 
field. We will consider two measurement schemes that can be used for this purpose: direct 
photon counting and random phase homodyne detection. 

The photon number distribution is given by the expectation value of a family of projec- 
tion operators p v = where \v) is the i^th Fock state. In principle, knowledge of this 
distribution enables us to evaluate any phase-independent observable related to the mea- 
sured field. A simple yet nontrivial observable, which we will use to point out difficulties 
with the completeness of the reconstructed information on the quantum state, is the parity 
operator: 

oo 

n = E(-i)»H (is) 

This operator is bounded, and well defined on the complete Hilbert space of a single light 
mode. Its expectation value is given by the alternating series of the photon number distri- 
bution: 

oo 

<n> = £(-!)"<&>, (16) 

which is absolutely convergent for any quantum state. Therefore, any pathologies connected 
to its determination from experimental data, if there appear any, cannot be ascribed to its 
singular analytical properties. 



A. Direct photon counting 

First, we will consider the reconstruction of phase insensitive properties of a single light 
mode from data measured using a realistic, imperfect photodetector. The positive operator- 
valued measure p n describing the probability of ejecting n photoelectrons from the detector 
is given by plj : 



Pn = ■ ; — exp(-7?a'a) :, (17) 



TV. 



where a is the annihilation operator of the light mode, and rj is the quantum efficiency of 
the photodetector. In the limit rj —>■ 1 we get directly p n = \n)(n\. In a general case, the 
probability distribution for the photoelectron number is related to the photon statistics via 
the Bernoulli transformation. This relation can be analytically inverted p^fl , which yields 
the expression: 
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P, = E^Pn, (18) 
ra=0 

where the kernel functions are given by: 

rSS = < ' M V IV'"" . . , (10) 




n > v. 



The inversion formula has a remarkable property that p v depends only on the "tail" of the 
photocount statistics for n > v. It has been shown that the inverse transformation can 
be applied to experimentally determined photocount statistics for an arbitrary state of the 
field, provided that the detection efficiency is higher that 50% P2| . 

We will now discuss statistical properties of the photon number distribution determined 
by photon counting within the general framework developed in Sec. [HI The case of perfect 
detection is trivial for statistical analysis. Therefore we will consider nonunit detection ef- 
ficiency, which is numerically compensated in the reconstruction process using the inverse 



Bernoulli transformation according to Eq. (18). For all examples presented here, the effi- 



ciency is r] = 80%, which is well above the 50% stability limit. 

In Fig. [I] we depict the reconstructed photon number distributions for a coherent state, 
a thermal state, and a squeezed vacuum state. The mean values E(p„) along with their sta- 
tistical errors [Var(p J/ )] 1//2 are compared with Monte Carlo realizations of a photon counting 
experiment, with the number of runs N = 4000. It is seen that for a thermal state and a 
squeezed vacuum state, the statistical error of the probabilities p v grows unlimitedly with 
the photon number v. However, any experimental histogram obtained from a finite number 
of runs ends up for a certain count number, and therefore the reconstructed photon statistics 
is zero above this number. An important feature that is evidently seen in the Monte Carlo 
simulations, are correlations between the consecutive matrix elements. The reconstructed 
photon number distribution clearly exhibits oscillations around the true values. This prop- 
erty can be quantified using the correlation coefficient defined in Eq. (0), which we plot for 
all three states in Fig. |^. For large u's, Corr(p I/ , p u +i) is close to its minimum allowed value 
— 1, which acknowledges that statistical correlations are indeed significant. 

These correlations affect any quantity computed from the reconstructed photon number 
distribution. The parity operator is here a good example: since in Eq. (|16D we sum up 
consecutive p„'s with opposite signs, their statistical deviations do not add randomly, but 
rather contribute with the same sign. Consequently, the statistical error of the evaluated 
parity operator may be huge. It is therefore interesting to study this case in detail. In prin- 
ciple, we could obtain statistical properties of the reconstructed parity using the covariance 
matrix for the photon distribution. However, it will be more instructive to express the parity 
directly in terms of the photocount statistics, and then to apply the statistical analysis to 
this reconstruction recipe. This route is completely equivalent to studying evaluation of the 
parity via the photon number distribution, as all transformations of the experimental data, 
which we consider here, are linear. 

A simple calculation combining Eqs. (|15| ) and ([18|) shows that: 




11 = E 1 " - Pn, (20) 
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where we have introduced in the upper summation limit a cut-off parameter K for the pho- 
tocount number. This formula clearly demonstrates pathologies related to the determination 
of the parity operator. For any 77 < 1, the factor (1 — 2/i]) n is not bounded, which makes 
the convergence of the whole series questionable in the limit K — > 00. Of course, for an 
experimental histogram the summation is always finite, but the exploding factor amplifies 
contribution from the "tail" of the histogram, where usually only few events are recorded, 
and consequently statistical errors are significant. 

Let us study these pathologies more closely using examples of a coherent state and a 
thermal state. For a coherent state \a), both the expressions for E(I1) and Var(Il) are 
convergent with K — > 00. However, the variance, given by the formula: 



Var(n coh ) = I 



/4(1 -77) . |2 \ . A . l2 ; 
exp — - a - exp(-4 a ) 

V V J 



(21) 



grows very rapidly with the coherent state amplitude a, when the number of runs N is 
fixed. For a thermal state with the average photon number n, the matter becomes more 



delicate. When K — > 00, the series (|20|) is convergent only for n < 1 / [2 ( 1 — 77)], which for 
1] = 80% gives just 2.5 photons. Even when the mean value exists, the variance is finite only 
for n < r\j [4(1 — rf)] and equals: 

= V ( ,-4-„) - (TtW) ' (22) 

We illustrate these results with Fig. || depicting Monte Carlo simulations for various average 
photon numbers. For coherent states, statistical fluctuations can in principle be suppressed 
by increasing the number of runs. For thermal states, the situation is worse: when n > 1, 
the variance cannot even be used as a measure of statistical uncertainty. 

Let us recall that the bound rj > 50% for the stability of the inverse Bernoulli transforma- 
tion is independent of the state to be measured. It has been obtained from the requirement 
that in the limit K — > 00 both E(p u ) and Var(p^) should converge. The example with 
the parity operator clearly shows, that the condition r\ > 50% does not guarantee that the 
reconstructed photon number distribution can be safely used to determine an arbitrary well- 
behaved phase independent observable. Thus, imperfect detection is inevitably connected 
with some loss of the information on the measured quantum state. 

We have noted that as long as finite, experimental data are concerned, evaluation of 
observables via intermediate quantities is equivalent to expressing them directly in terms of 
the measured probability distributions. One might try to circumvent the 77 > 50% bound for 
reconstructing the photon statistics by applying the inverse Bernoulli transformation in two 
or more steps, and compensating in each step only a fraction of its inefficiency Of course, 
such a strategy must fail, as for any finite sample of experimental data this treatment 
is equivalent to a single transformation which is unstable. In many-step processing this 
instability would be reflected in increasing correlations and statistical errors exploding to 
infinity. 
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B. Random phase homodyne detection 



Random phase homodyne detection is a recently developed technique for measuring 
phase-independent properties of optical radiation, which goes beyond certain limitations 
of plain photon counting ||23[j . Data recorded in this scheme is the difference of counts on 
two photodetectors measuring superposition of the signal field with a strong coherent local 
oscillator. The count difference is rescaled by the local oscillator amplitude, and the resulting 
stochastic variable x can be treated with a good approximation as a continuous one. The 
photon number distribution is reconstructed from the random phase homodyne statistics 
p(x) by integrating it with pattern functions f v (x) 

(p u }= dxf u (x)p(x). (23) 



A convenient method for numerical evaluation of the pattern functions has been described 
in Ref. §. 

Let us now discuss statistical properties of the homodyne scheme in its discretized version 
used in experiments, when the rescaled count difference is divided into finite width bins. As 
the local oscillator phase is random, the setup has no controllable parameters, and the 
statistics of the observables is fully determined by p(x). Statistical errors of the density 
matrix in the Fock basis reconstructed via homodyne detection have been studied in Ref. || . 
Here we will focus our attention on statistical correlations exhibited by the diagonal density 
matrix elements, and their further utilization for evaluating phase-independent observables. 

We will consider the unit detection efficiency r] = 1, with no compensation in the pro- 
cessing of the experimental data. This is the most regular case from the numerical point of 
view. When 77 < 1 and the compensation is employed, the statistical errors are known to 
increase dramatically M. In Fig. |] we depict the homodyne reconstruction of the photon 
number distribution for the three states discussed in the previous subsection. For large z/, 



the statistical errors tend to a fixed value y2/N, which has been explained by D'Ariano et 
al. using the asymptotic form of the pattern functions ||. Again, Monte Carlo simulations 
suggest that the reconstructed density matrix elements are correlated, which is confirmed by 
the correlation coefficient for the consecutive photon number probabilities, plotted in Fig. ||. 
A simple analytical calculation involving the asymptotic form of the pattern functions shows, 
that for large v this coefficient tends to its minimum value —1. 

One may now expect that no subtleties can be hidden in using the reconstructed photon 
number distribution to evaluate the parity operator according to Eq. ([16]). However, let 
us recall that the parity operator is equal, up to a multiplicative constant, to the Wigner 
function at the phase space origin |H|j25| . The Wigner function is related to the homodyne 
statistics via the inverse Radon transformation, which is singular. In particular, applying 
this transformation for the phase space point (0, 0) we obtain the following expression for 
the parity operator in terms of the homodyne statistics 



(fl) = 1 dxp(x)^-P-, (24) 

Z J —00 CLX X 

where P denotes the principal value. Due to the singularity of the inverse Radon transform, 
its application to experimental data has to be preceded by a special filtering procedure. This 
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feature must somehow show up, when we evaluate the parity operator from the reconstructed 
photon statistics. In order to analyse this problem in detail let us discuss evaluation of the 
truncated parity operator tlx from a finite part of the photon number distribution: 

K roo 

(IM = £(-ir<A,> = / g K (x)p(x), (25) 

where 

9k{x) = J2(-irUx) (26) 

can be considered to be a regularized kernel function for the parity operator. In Fig. |6| we plot 
this function for increasing values of the cut-off parameter K. It is seen that the singularity 
of the kernel function in the limit K — ► oo is reflected by an oscillatory behaviour around 
x = with growing both the amplitude and the frequency. This amplifies the statistical 
uncertainty of the experimental homodyne data. In Fig. [7| we show determination of the 
parity operator for the three states discussed before, using increasing values of the cut-off 
parameter K. Though we are in the region where the true photon number distribution 
is negligibly small, addition of subsequent matrix elements increases the statistical error 
in an approximately linear manner. This is easily understood, if we look again at the 
reconstructed photon number distributions: increasing K by one means a contribution of the 
order of J 2/N added to the statistical uncertainty, and, moreover, these contributions tend 
to have the same sign due to correlations between the consecutive matrix elements. Thus, 
determination of the parity operator from homodyne statistics requires an application of a 
certain regularization procedure. It may be either the filtering used in tomographic back- 
projection algorithms, or the cut-off of the photon number distribution. The statistical 
uncertainty of the final outcome is eventually a result of an interplay between the number 
of experimental runs and the applied regularization scheme. 

Finally, let us briefly comment on the compensation for the nonunit efficiency of the 
homodyne detector. First, one might think of applying a two-mode inverse Bernoulli trans- 
formation directly to the joint count statistics on the detectors. However, it is impossible 
in the homodyne scheme to resolve contributions from single absorbed photons due to high 
intensity of the detected fields. The inverse Bernoulli transformation has no continuous 
limit, as consecutive count probabilities are added with opposite signs. Nevertheless, the 
nonunit detection efficiency can be taken into account in the pattern functions ||. In this 
case the statistical errors increase dramatically, and explode with v — > oo, which makes de- 
termination of the parity operator even more problematic. This is easily understood within 
the phase space picture: the distributions measured by an imperfect homodyne detector are 
smeared-out by a convolution with a Gaussian function P7PS| . Evaluation of the parity 
operator, or equivalently, the Wigner function at the phase space origin requires application 
of a deconvolution procedure, which enormously amplifies the statistical error |29|| . 

IV. CONCLUSIONS 

We have presented a complete statistical analysis of determining quantum observables 
in optical measurement schemes based on photodetection. We have derived an exact ex- 
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pression for the generating function characterizing statistical moments of the reconstructed 
observables, which, in particular, provides formulae for statistical errors and correlations 
between the determined quantities. These general results have been applied to the detection 
of phase-independent properties of a single light mode using two schemes: direct photon 
counting, and random phase homodyne detection. This study has revealed difficulties re- 
lated to the completeness of the reconstructed information on the quantum state: in some 
cases the parity observable, which is a well-behaved bounded operator, effectively cannot 
be evaluated from the reconstructed data due to the exploding statistical error. 

We have recalled that the parity operator is directly related to the value of the Wigner 
function at the phase space origin. Thus, our example can also be interpreted as a particular 
case of the transformation between two representations of the quantum state: in fact, we 
have considered evaluation of the Wigner function at a specific point (0, 0) from the relevant 
elements of the density matrix in the Fock basis. Therefore, our discussion exemplifies 
subtleties related to the transition between various quantum state representations, when 
we deal with data reconstructed in experiments. Though a certain representation can be 
determined with the statistical uncertainty which seems to be reasonably small, it effectively 
cannot be converted to another one due to accumulating statistical errors. 

The presented study suggests that the notion of completeness in quantum state mea- 
surements should inherently take into account the statistical uncertainty. From a theoret- 
ical point of view, the quantum state can be characterized in many different ways which 
are equivalent as long as expectation values of quantum operators are known with perfect 
accuracy. In a real experiment, however, we always have to keep in mind the specific experi- 
mental scheme used to perform the measurement. This scheme defines statistical properties 
of the reconstructed quantities, and may effectively limit the available information on the 
quantum state. Determination of a family of observables does not automatically guarantee 
the feasibility of reconstructing the expectation value of an arbitrary well behaved operator. 
Reconstruction of any observable should be preceded by an analysis how significantly the 
final result is affected by statistical noise corrupting the raw experimental data. 
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FIG. 1. Reconstruction of the photon number distribution from photon counting for (a) a 
coherent state with (h) = 4, (b) a squeezed vacuum state with (n) = 1, and (c) a thermal state 
with (n) = 2, from N = 4000 runs in each case. Monte Carlo simulations of a photon counting 
experiment, depicted with points, are compared with exact values (solid lines), with the statistical 
errors [Var(/o jy )] 1 / 2 marked as grey areas. The detection efficiency is r] = 
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FIG. 2. The correlation coefficient between the consecutive density matrix elements 
Corr(/9 l/ , Pv+i), depicted for the coherent state (o), the squeezed state (□), and the thermal state 
(*) from Fig. ffl. 




FIG. 3. Determination of the parity operator for (a) coherent states and (b) thermal states 
with the increasing average photon number. Each square represents the parity evaluated from 
Monte Carlo simulated photon statistics with N = 4000 runs for a given average photon number. 
The solid lines and the grey areas depict the mean value E(n) and the error [Var(II)] 1 / 2 . 
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FIG. 4. Random phase homodyne reconstruction of the photon number distribution for (a) the 
coherent state (b) the squeezed state and (c) the thermal state, all states with the same photon 
numbers as in Fig. g. The range of the homodyne variable is restricted to the interval — 6 < x < 6 
divided into 1200 bins. The simulated homodyne statistics is obtained from N = 4 • 10 4 Monte 
Carlo events. 
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FIG. 5. The correlation coefficient for the consecutive photon number probabilities for the 
coherent state (o), the squeezed state (□) and the thermal state (*). 




FIG. 6. Regularized kernel functions for the parity operator gxix) evaluated as a finite sum of 
the Fock states pattern functions, for increasing values of the cut-off parameter K. 
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FIG. 7. Reconstruction of the truncated parity operator Tlx for the coherent state (o), the 
squeezed state (□), and the thermal state (*) with various values of the cut-off parameter K, 
using the same Monte Carlo homodyne statistics as in Fig. The simulations are compared with 
the corresponding mean values E(IIjf) and errors [Var(IlK-)] 1 ' 2 plotted as solid lines surrounded 
by grey areas. 
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